import numpy as np
import matplotlib.pyplot as plt
from sklearn.preprocessing import StandardScaler

# ⑧　完成数据的读取、初始化和特征缩放等相关操作    （4分）
# load
d_train = np.loadtxt(r'../../../../../large_data/finalX/train_chem.txt', delimiter=',')
m_train = len(d_train)
d_test = np.loadtxt(r'../../../../../large_data/finalX/test_chem.txt', delimiter=',')
m_test = len(d_test)
# scale
scaler = StandardScaler()
d_train = scaler.fit_transform(d_train)
d_test = scaler.fit_transform(d_test)
# init
X_train = np.c_[np.ones(m_train), d_train[:, :-1]]
y_train = d_train[:, -1]
X_test = np.c_[np.ones(m_test), d_test[:, :-1]]
y_test = d_test[:, -1]


# ⑨　正确自定义costFunction，并正确传递函数参数     （4分）
def costFunction(h, y):
    m = len(h)
    e = h - y
    return 1 / 2 / m * e.T.dot(e)


# ⑩　实现梯度下降函数，每一次迭代输出当前的代价函数值   （4分）
def grad(X, y, alpha=0.01, iter0=2000):
    m, n = X.shape
    group = iter0 // 20
    theta = np.zeros(n)
    j_his = np.zeros(iter0)  # cost function value history in iterations
    for i in range(iter0):
        h = X.dot(theta)
        j = costFunction(h, y)
        j_his[i] = j
        if 0 == i % group:
            print(f'#{i + 1} cost function value = {j}')
        dt = 1 / m * X.T.dot(h - y)
        theta -= alpha * dt
    if 0 != i % group:
        print(f'#{i + 1} cost function value = {j}')
    return theta, j_his, h


def score(h, y):  # ATTENTION
    u = np.sum((h - y) ** 2)
    mu = y.mean()
    v = np.sum((y - mu) ** 2)
    return 1 - u / v


# 11　执行梯度下降函数，并使用matplotlib画出代价随迭代次数变化的趋势    (4分)
plt.figure(figsize=[16, 8])
spr = 1  # subplot row
spc = 2  # subplot column
spn = 1  # subplot number
plt.subplot(spr, spc, spn)
alpha=0.01
iter0=2000
theta, j_his, h_train = grad(X_train, y_train, alpha, iter0)
print(f'Theta = {theta}')
print(f'Training score = {score(h_train, y_train)}')
plt.plot(j_his, label='cost function')
plt.grid()
plt.legend()
plt.xlabel('Iterations')

# 12　计算模型在test_chem.txt的评估指标，并输出     4分
h_test = X_test.dot(theta)
score_test = score(h_test, y_test)
print(f'Testing score = {score_test}')

# 13　使用matplotlib画图：x轴为真实值，y轴为预测值，画出测试集的真实值和预测值对比图，并在标题中标出精度        (4分)
spn += 1
plt.subplot(spr, spc, spn)
plt.scatter(y_test, y_test, s=1, zorder=0, label='target value')
plt.scatter(y_test, h_test, s=1, zorder=10, label='hypothesis value')
plt.title(f'Score = {score_test:.2f}')
plt.grid()
plt.legend()
plt.xlabel('real value')

# show all drawings
plt.show()
